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Abstract 

We demonstrate an iterative scheme to approximate the optimal trans- 
portation problem with a discrete target measure under certain standard 
conditions on the cost function. Additionally, we give a finite upper bound 
on the number of iterations necessary for the scheme to terminate, in terms 
of the error tolerance and number of points in the support of the discrete 
target measure. 

1 Introduction 

In p] , Caffarelli, Kochengin, and Oliker developed a numerical algorithm to cal- 
culate approximations to discretizations of the far-field reflector antenna prob- 
lem, and also gave a finite upper bound on the number of steps necessary. Later, 
Xu-Jia Wang showed that the far-field reflector problem can be interpreted as 
an optimal transportation problem (see |10j . and also [6]). A number of similar 
algorithms have also been considered, each for a problem related to some form of 
optimal transportation problem: [5] proposes an algorithm to solve the classical 
Monge- Ampere equation in two dimensions (see however, Remark 13.21 below), 
while [5] considers the near-field reflector problem where the target domain is 
contained in a flat domain. Such results lead to a natural question: can a similar 
numerical algorithm can be developed for more general optimal transportation 
problems (and perhaps more importantly, with a similar upper bound on the 
number of iterations)? 

In this paper, we answer this question in the affirmative, given certain con- 
ditions on the optimal transportation cost function. Namely, we show that an 
iterative scheme similar to the one given by Caffarelli, Kochengin, and Oliker 
can be applied to optimal transportation problems associated to a cost function 
satisfying the conditions ( |Reg[ ), (|Twist[) . ( |Nondcg[ ), and (MTW + ) (one can ver- 



ify that the optimal transportation formulation of the far-field reflector problem 
satisfies the required conditions, see Section [5] for relevant definitions). Addi- 
tionally, we give a finite upper bound on the number of iterations necessary 
for this scheme to halt, which is of the same form as what is proven in pQ. 
We are careful to note, however, that we do not provide a bound on the actual 
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computational complexity of the scheme, which may vary depending on the geo- 
metric details of the particular cost function (towards this direction, Kochengin 
and Olikcr further explore the actual implementation of the algorithm for the 
far-field reflector problem in [J]). 

The aforementioned conditions are natural in the theory of optimal trans- 
portation, (MTW+I was first introduced by Ma, Trudinger, and Wang in [7J, 
where they show it is a key condition leading to a proof of regularity properties. 



Later, Loeper proved in [5] that a weaker form of (MTW+ ) is actually necessary 



for regularity, and also showed a number of geometric implications of ( MTW 



and its weaker form (see also [2]). It is precisely these geometric implications 
that we exploit in this paper, namely that certain sublevel sets possess a gener- 
alized notion of strong convexity, and use this to show that calculations similar 
to Caffarelli, Kochengin, and Oliker can be applied to this more general case. 

The organization of this paper is as follows. In Section [5] wc introduce the 
optimal transportation problem, and a number of basic conditions and results 
classical to the literature. In Scction|3]we state the main results that are claimed 
in this paper, while in Section 0] we give a detailed description of the aforeme- 
tioned iterative scheme. Section [5] is devoted to a series of calculations, which 
are invoked in Section [5] in order to provide the claimed upper bound (|3.1I) on 
the number of iterations necessary. 



2 Preliminaries of Optimal Transportation 

In this section, we recall some basic facts and definitions regarding the optimal 
transportation problem, along with the key conditions introduced in [TJ. For a 
(much) more comprehensive study of the optimal transportation problem, the 
interested reader is referred to [pj. 

Given two domains fl and Cl, and two probability measures fi and v defined 
on them, along with a real valued cost function c defined on f2 cl x Q , we wish 
to find a measurable mapping T : — > 51 satisfying T#[i = v (defined by 
T#(i (E) = fj, (T _1 (.E)) for all measurable IcO) such that 

/ c(x,T(x))dfi — min / c(x, S(x))dfi. 
Jn s #v=vJn 

We will say that such a T is a solution to the optimal transportation problem, 
or a Monge solution, for the cost function c transporting the measure // to the 
measure v. 

Under mild conditions on c and the measure [i, it is known that a unique 
solution to this problem exists. For instance, it is sufficient if Q and Cl are 
bounded subsets of Riemannian manifolds M and M respectively, c satisfies 
conditions ( |Rcg[ ) , (|Twist[) . and ( |Nondcg| below, and fi is absolutely continuous 
with respect to the volume measure dVolju defined by the Riemannian metric 
on M. Additionally, in this case T can be determined from a scalar valued, 
Lipschitz continuous potential function (f> by the formula 

T(x) = expl(D(f>(x)) 
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where the differential D<f>(x) exists in the dVolM almost-everywhere sense, and 
exp x (p) is defined below (see [H Chapter 10]). 

Now let fl and Si be open, bounded domains with piecewise smooth bound- 
aries in (i-dimensional Riemannian manifolds M and M respectively. We intro- 
duce a number of key conditions on the cost function c. Below, Dc and Dc are 
the differential of c in the x and x variable respectively, which unambiguously 
defines the mapping 

-DDc(x,x) : T x n -> T^ c[x ^ (T|Sl) = T*S1. 



Smoothness of cost function: 

Assume that 

ce(7 4 (!l cl xO cl ). (Reg) 

Twist condition: 

c satisfies ([Twist [I if for each xq £ Sl cI and xq £ S7 cl , the mappings x H> 
— Dc{xo,x) and x i-> — Dc(x, xo) are injective. 

For any p £ —Dc(xq, fi) and xo S S7, (resp. p S — -Dc(S7, So) and So 6 S7) we 
write exp Xo {p) (resp. exp% Q {p)) for the unique element of S7 (resp. S7) such that 

-Dc(x ,exp c {p)) =p, 

n / c / - n (Twist) 
--Dc(exp So (p),x ) = p. 

Remark 2.1. For fixed x 6 S7 and 5 6 S7, we will denote the representations of 
sets E C S7 and £7 C S7 in the cotangent spaces above x and x by 

: = -£<£,*), 
[£] i : = -flc(:c 1 £). 

Nondegeneracy condition: 

c satisfies (Nondeg) if the following linear mapping is invertible for every 
x £ ST 1 and x £ ST 1 : 

- DDc(x, x) : T x fl -> T|fL (Nondeg) 
Strong MTW condition: 



A cost c satisfies (MTW+J) if there exists some So > such that 

- (c*„m - Ci jtPt r''c^' k c«' l {x, x)V^ mm > <5o|y| 2 |?7| 2 , (MTW+) 

for any x £ ST 1 , x £ ST 1 , and V £ T x fl and 77 e T*S1 such that rj(V) = 0. Here 
all derivatives are with respect to a fixed coordinate system, regular indices 
denote derivatives of c with respect to the first variable, while indices with a 
bar above denote derivatives with respect to the second derivative, and a pair 
of raised indices denotes the matrix inverse. 

Additionally, we need the following concepts of c-convexity of a domain, and 
c-convex functions. 
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Definition 2.2. We say that a set E C fl is (strongly) c-convex with respect to 
xo G fi if the set is a (strongly) convex subset of T? o Q. 

Similarly, we say E C f2 is (strongly) c-convex with respect to xq G fi if the 
set \E\ x is a (strongly) convex subset of T* Q tt. 

Remark 2.3. For some fixed x E fi, given any two points p\ := —Dc(xi,x) and 
p 2 := — Dc(x2,x) with x±, x 2 G fi, we define the c-segment with respect to x 
between x\ and x 2 as the image of the straight line segment between pi and pi 
under the map exp%{-). It is clear that E is c-convex with respect to x G fi if 
and only if every c-segment with respect to x between any two x\ and xi G E 
remains inside E. An symmetric statement and definition holds with the roles 
of fi and fi reversed. 

Definition 2.4. We say that a function <f> is c-convex, if for every xq G ft there 
exist xo £ Ct and Ao £ K such that 

<p(x Q ) = -c(x ,xo) + A , 

(p{x) > -c(x,X ) + Xq 

for all x xq. 

We also call such a function — c(-, xq) + Xo that satisfies the above equality 
and inequality, a c- support function to (f> at Xq. 



3 Main Result 

We state in this section, the main result of the paper. 

Assume that fi := IdVoiM for some positive, real valued / £ C°°(f2 cl ) which 
satisfies 



I(x) dVol M (z) = l, 



while the domains fi and f2 are bounded subsets of Riemannian manifolds M 
and M respectively. Additionally, assume the cost function c satisfies con- 



ditions peg]), (jTwistp . (|Nondcg[), and (|MTW + [), fix an integer K > 2 and 



Cj i7, a finite collection of K distinct points, and {fi}fLi, a collection of 



K real numbers satisfying ^2f =1 fi = 1 and < fi < 1. At this point, we make 
the additional assumption that ft is c-convex with respect to {xi}fL 1 and f2 is 
c-convex with respect to f2 (however, we do not make any assumptions on the 
support of /). 

The main result we present here is as follows: for any e > 0, there is an 
iterative scheme to find numbers /i > and d i G M such that 

\fi-fi\<s, Vl<i<K, (3.1) 
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and T{x) = exp%{D4>{x)) is the Monge solution for the cost function c trans- 



it 



porting the measure /i to the measure fi^x t , where 



max 

Ki<K 



-c(x, Xi) + d i 



(3.2) 



Additionally, we show that the number of iterations n e necessary to find {/f}|£i 
and {d i }fL 1 is bounded above: 



n, < K 



KCM sup J2 1 
SA 



where 



C : 



max sup 

i<**<* ae n \(-DDc(x,Xi)) (-Dc(x,Xi) + Dc{x,x k ))\ T ,_a 



max V. d ~ 1 (d[QL ) 

l<i<K v Xl ' 

det (-D5c(a;,Xi)) 



(3.3) 



M:= max sup e -c(*,s fc )+<<x,xi) + 1; 

Kk<K xen 

A:= min inf e -c(*,s fc )+c(x^O 



(5 : = 



K-V K 



(3.4) 



and H d 1 (■) denotes (d — 1) -dimensional Hausdorff measure. 
Remark 3.1. The upper bound here has the same dependency on the number 
of points K in the target and the error e as the one given in Q]. Note that 
the bound we give here is on the number of iterations (the number of c-convex 
functions <p n necessary, see the description in Section 2] below) . However, since 
there are exactly K intermediate steps per iteration, counting the total number 
of steps we would multiply by an additional factor of K. By the choice of 5, we 
note that the total number of steps with accuracy e has an upper bound that 
is essentially of the order of K A /e. 

Remark 3.2. We would like to point out here that it is possible to obtain some 



finite upper bound on n e under the degenerate version of ( MTW+ ) where the 
constant 5 Q 

|2 



= (this is the case, for example, with the cost function c(x, x) = 
\x — y\" on ill™, which is known to correspond to the classical Monge- Ampere 
equation). However, the difference will be that the higher order terms in the 
expression of G'(d + Ad) — G*(d) may only have order 0(|Ad|) instead of 
o(|Ad|) as | Ad | — > (see the proofs of Propositions 15.51 and 15 . As a result, 
the corresponding terms in the derivative of G l may not vanish, but can be 
bounded above, and the final upper bound may contain Y)^-i (^) more terms. 



4 Description of the iterative scheme 

In describing the algorithm, we follow most of the notation used by Caffarelli, 
Kochengin, and Oliker in Q] to better highlight the parallels between the two 



5 



schemes. 

Definition 4.1. Let <f> be a c-convex function on f2. Then for any measurable 
set E C fi, we define the visibility set of ui associated to <fi by 

V<f,(E) : = {x £ ^ | there exist A G K and y <E E such that, 
— c(-, x) + A is a c-support function to at x} . 

We also define 

G{4>,E) :=n{V^E)). 

Definition 4.2. Write d := (d 1 , d 2 , . . . , d K ) for d 1 > 0, 1 < i < K (we will use 
the shorthand d > to denote this from now on) . Then we define the following 
function: 

4>d{x) := max [— c(x, Xj) — log cf]. 

l<i<K 

We will generally use superscripts to denote coordinates for the remainder 
of the paper. 

Since we are only concerned with the optimal transportation problem when 
v is a finite sum of delta measures, we make some notational simplifications as 
follows. For any 1 < i < K we will write 

Vd(xi) := V^{{xi}), 
G i {4>) :=G{<f>,{xi}), 
G*(d) :=G(0 d ,{xJ), 

for any c-convex function <f> and d > 0, when the collection of points {xi}^ =1 is 
fixed. 

Remark 4.3. Since c satisfies conditions ( |Rcg[ ) and (|Twist[) , <j>d is differcntiablc 
dVolA/-a.e on f2, and we can define the map 

T d {x) := exp c x {D<t> d (x)) 

for dVolM-a.e. x e fi. Additionally, we can see that 

T d (x) = Xi 

for all x G Vd(xi), and since 4> d is clearly c-convex, this implies that T d is 
the unique Monge solution from the measure JdVolj\f to J2i=i G l {d)5 Si ([9, 
Chapter 10]). Intuitively the quantities G*(d) give the amount of mass in f2 
that is transported by the map T d to the point Xi for each 1 < i < K. 

We will make also use of one of the consequences of a theorem proved by 
Locper (0 Theorem 3.7]). 

Theorem 4.4 (Loeper). If c satisfies ( |Reg| , (jTwistj) . ( |Nondeg[ ), and ( |MTW + [) ; 

and f2 is c-convex with respect to f2, then for any c-convex function <f> and jg!!, 
Vrf,{x) is c-convex with respect to x. 
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Before describing the iterative scheme, we need to prove a monotonicity 
property of the functions G % . To do so, we first give an alternate characterization 
of the sets V d {x,i). 

Lemma 4.5. Suppose c, fi, £1, and Q satisfy all of the conditions of Section^ 
then for any 1 < i < K 

V d (xi) = V d>i 

if Va,i 7^ 0, where 

V d>i :={xen \ -c(x,Xi) -logd 1 > -c(x,x k )-logd k , VI < k < K) . 
Otherwise, 

ti(V d (xi)) = 0. 

Proof. Note that for any choice of d, <j) d is a c-convex function on S7. It is clear 
from the definition of (f> d that, at any point in V di i there exists a c-support 
function of the form — c(-, Xi) + A supporting to <p d at that point, thus we obtain 

V d ,i C V d (xi). 

If V d (xi) = 0, the claims of the lemma immediately follow, so let us suppose 
there exists some x$ € V d (xi). Then there exists some A G R such that 

-c(x Q ,Xi) + A = (f> d (x ) 
-c(x,Xi) + A < 4> d (x) 

for all x £ CI. In particular by the definition of 4>d, 

-c(x ,Xi) + A = (j> d (x ) > -c(x ,Xi) - logd 1 
A > - logd 1 . 

If there exists a point x\ £ V d ^, we have from the definition of <p d and V d ,i 
that 

-c(xi,Xi) + A < (/) d (xi) = -c(xi,Xi) - logd 1 
=^ A < -logtf, 

hence A = — log d l and we obtain xq £ V dl i- 

Now suppose Vd,i = 0- Assume first that 4>d is differentiable at Xq. By the 
definition of </>d, this implies that for some 1 < jo < k, i ^ jo we have 

4>d(x ) = -c(x ,x jo ) - logd 30 , 
(j)d(x) > -c(x,x JO ) -logd 10 , V x ^ x 

and in particular 

D4> d (xo) = -Dc(x ,x ]Q ). 
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At the same time, since xq £ Vd{xi), we must have 

D(j> A {x Q ) = -Dc(x Q ,Xi). 

which contradicts (|Twist|) . However, since <pd is diffcrentiable dVolM almost 
everywhere (Remark l4.3[) , and /x is absolutely continuous with respect to dVol^/, 
we obtain \x (Vd(xi)) = as desired. □ 

Now we can easily see the following corollary. 

Corollary 4.6. Suppose c, fi, Q, and f2 satisfy all of the conditions of Section^ 
Also, fix an index 1 < i < K , another index j ^ i, and values d k > for all 
k i. Then, G*(d) is decreasing in d? , G J (d) is increasing in d l , and we have 
the following limits: 

lim G l (d) = 1, (4.1) 
lim G^fd) =0. (4.2) 

d'-i-O 



Proof. If all d k are fixed for k ^ i, we can easily see that the sets Vd.i are 
increasing as d l decreases. Hence by Lemma 14.51 above we immediately obtain 
the following monotonicity property for d\ < d\: 

I 1 (v d d\ d* (■''•') > /' d\ d* (•'•■;) • 

Similarly, since j ^ i, we obtain when d\ < d\: 

I' fail,..., d*,.. < /' (V d i ,/, ,/<. i) • 

Then recalling that G'(d) = (i(Vd(xi)), the claims of monotonicity arc imme- 
diate. 

To obtain the limiting values, note that from the boundedness of c and the 
definition of <j>&, if all d k are kept fixed for k ^ i, we will have Vd,i = while 
Vd,j = for d l > sufficiently small. Then, again by Lemma 14.51 above, the 
limits (|4.1[) and (|4.2[) are immediate. □ 

With this monotonicity property in hand, we are ready to describe the iter- 
ative scheme in detail. 

If K = 1, we have that — c(-,xi) + d gives rise to an optimal solution from 
\x to for any choice of d G R (the associated mapping is simply T(x) = Si). 
If K = 2, we let 

^>(a;) := max{— c(x, Xi), —c{x,X2) + d} 

then adjust d continuously until /j := G 1 (d) and /a := G 2 (d) satisfy the desired 
bounds. Thus we assume that K > 3. 
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The algorithm now consists of starting with a c-convex function <fio := c/)^o 
for an appropriate choice of d°, then decreasing each of the parameters d k in 
turn until the desired mapping is found. 

First, (after fixing e > 0) we define 5 > by 



/ 



K - 1 ' K 

(the second restriction will play a role in showing the upper bound on the 
number of steps necessary, see Section [6]). Then, define the set $5 by 

$5 := {cf>d I d > 0, G l (d) < ft + 5 for all 2 < i < K) . 

Since fa + S > for each 1 < i < K, by the limiting value (|4.2[) in Corollary |4.6l 
we see that s $5 whenever d 1 > is sufficiently small, in particular $5 7^ 0. 
Take any element of §5 and let it be 4>o- We now construct a sequence of 
c-convex functions <p n 6 $5 as follows. 

Suppose we have n £ $5, we construct 0„+i € $5 by first constructing a 
sequence of K intermediate c-convex functions. Let 4> n ,i '■= 4>n- Then, for any 
1 < i < K — 1, suppose 



id 1 .,d l 



G *«< 



If I G l+ 1 (0TL,t ) — .A+ 1 < £ we simply set (j> n ,i+i ■= 4>n.i- Otherwise, since cf) n ,i <E $s 
and i + 1 > 2, we must have 

G' 1+1 (0, M ) < fi+i-S. 

Now, G t+1 (d) is continuous in d z+1 by Proposition 15 .61 in Section[5l and has the 
monotonicity property described in Corollary 14.61 Since /j+x < 1, and we have 
the limit (|4.1|) , we can find a < d^l < d l +] such that 

C t+1 (^n,iJ dn,i-> ■ ■ • i 'Cji J - J ^n,i) G + <*>)• 



Taking n) »+i := 4>^ d i d 2 gt+i d K y we see by Corollary 14.61 again that 
4>n,i+i G $5. We continue in this manner for each 1 < i < K — 1 until we 
determine 4> n +x '■= 4>n,K- 

If it happens that on the n E iteration we have 

4>n e = 4>n e ,l = ... = 4>n E ,K = 4>n c +l, 

this would imply that 

|G i (^ E )-/. i |<<5<-^ T <e 
for all 2 < i < K by the choice of 5. At the same time, 

K K 

IG 1 ^) - h\ = |1 - J2 G*(0„J - 1 + M 

i=2 i=2 

< (K - 1)6 < e. 



Thus we can see that if we let 

d t := - log 4^ 

the construction (J3T2J) gives us exactly the c-convex function cf> ne , while the fi 
satisfy the desired bound (|3.1|> . and Remark |4~B1 implies that </f>„ e gives rise to 
the Monge solution with the desired properties. Notice that the value of d 1 
remains fixed throughout the algorithm. 

Now, it is a priori possible that this scheme may continue for an infinite 
number of iterations. The following two sections are devoted to showing that 
this is not the case, and showing the upper bound (|3.3[) on the number of 
iterations. 

5 Derivative of the map G 

We now show that G l is differentiable in the ith variable, and obtain an upper 
bound for this partial derivative. This bound will be crucial in showing the 
desired upper bound (|3.3p . 

Throughout this section, we will fix d > and one particular index 1 < i < 
K. Now define 

Ad := (Ad 1 ,..., Ad*) e M A ', 
|Ad| := max \Ad k \, 

\<k<K 

with Ad 1 = 0. It will be implicitly assumed that d k + Ad k > for each 
1 < k < K. 

We will also write for any 1 < j < K, 

Wi,i,d = Wj := {x e | — c(x,Xj) - logd 7 < — c(x,Xi) - logd 1 } , 
W££ d = d :={xeil \ -c(x, xj) - log {d 1 + Ad 1 ) < -c(x, x t ) - log d 1 } 

Remark 5.1. Note that Wj is the visibility set V<f,(xi) for the c-convex function 
<f)(x) := sup {— c(x, x~i) — logd 1 , — c(x, Xj) — logd J }. Hence, by Theorem 14.41 we 
see that for any 1 < j < K, Wj is c-convex with respect to X{. 

We define these sets here because when we consider various difference quo- 
tients of G\ we will obtain intersections of sets of the form W^ A \ Wj with 
different indices. 

First, a technical lemma. The author believes this is a well known fact, but 
in the interest of completeness a proof is provided here. 

Lemma 5.2. If A C B are both bounded, convex sets in M. d , then 

W 1 - 1 (OA) < U d - 1 (dB) 
where 'H d ~ 1 (•) is (d — 1) -dimensional Hausdorff measure. 
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Proof. First, if the affine dimension of A is strictly less than d—1 then (dA) = 
and the claim is immediate. If the afhnc dimension of A is d, since A is convex, 
for each p £ dA\N there is a unique unit vector v(p) such that (p' —p,v(p)) < 
for all p' <E A, where iV C dA satisfies (iV) = 0. If the affine dimen- 

sion of A is d — 1, we fix one of the unit vectors that is normal to the d—1 
dimensional affine hull of A, and choose v{p) to be that vector (in this situa- 
tion, N = 0). In both cases, since A C B and B is convex and compact, if 
we define \{p) := sup {A > | p + Xv(p) £ B}, we see that X(p) is finite and 
p + X(p)v(p) e dB. Thus, 

*(p) A(p)u(p) 

is a well-defined map from dA \ N to 95. 

We now claim that * is injective on dA \ N, and also 

1*00 - ^(p 2 )| > bi -p 2 | (5.1) 

for any pi, p^ G cM \ N. Since (A^) = 0, by the definition of Hausdorff 

measure we see that this would prove the lemma. 

To show injectivity, suppose p\ ^ p2 are both in dA \ N but 

= *(p 2 ) =: go- 
Then, 

\pi - <7o I 2 = \P2 - *(pi)| 2 

= \P2 - (Pi + X( Pl )v( Pl ))\ 2 

= \P2 -Pi\ 2 - 2A(pi)(p 2 ~Pi,v(pi)) + A(pi) 2 
> X( Pl ) 2 

= bi-*(pi)| 2 
= \Pi - <7o| 2 - 

However, by reversing the roles of p\ and P2 above we obtain the opposite strict 
inequality, hence VP must be injective. 

Now to prove the expansion property, first note that for any pi, P2 G dA\N, 

(P2 -pi,*(pi) - Pi) = X(pi)(p2 ~Pl,v(pi)) 

< o 

and similarly 

(pi - P2, ^(P2) - P2) < 0. 
By adding these two equations together, we obtain 

> (P2 -Pl,*(Pl) ~Pl) + (P2 -Pl,P2 ~ *(P2)} 
= (p 2 -Pl,*(Pl) - *(P2)) + \P2 - Pi | 2 , 
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and by applying the Cauchy-Schwarz inequality, 



|P2 - Pi| 2 < {p2-Pi,*(Pa)-*(pi)) 
< b2-j»i||*(pa)-*(pi)| 

which readily implies inequality (|5 . 1 [) . □ 



We now prove a key lemma, utilizing the (MTW+ I condition. Here, we are 



exploiting the strong convexity of sublevelsets of c-functions under ( |MTW_ 



which was proven originally by Loeper. Below, djvf (•, •) and d M (-,-) are the 
geodesic distances given by the Riemannian metrics on M and M respectively. 

Lemma 5.3. Suppose [i, c, ; and CI satisfy all of the conditions in Section^ 
If xq £ f2 cl is a point where 

-C(XO) Xi) + c(x , Xj) = -c(x , Xi) + c(x ,x k ) 

—Dc(xo, Xi) + Dc(x , Xj) — a (~Dc(xo, Xi) + Dc(x , Xk)) 

for some a ^ and j ^ k, then there exists some r > depending on i, j , and 
k such that at least one of either 

-c(x,Xi) + c(x, Xj) ^ —c(xo, Xi) + c(x ,Xj) 

or 

-c(x, Xi) + c(x, x k ) ^ -c(x , Xi) + c(x , x k ) 
is true for all < (1m (x, xq) < r. 

Proof. By condition (|Twistj) . we find that ct^l. Thus, (writings := — Dc(xq, Xi), 
Pj := — Dc(xq, Xj), andpfc := —Dc(xq, Xk), which are all distinct again by (jTwistl ) 
we see that pi, pj, and p k are all collinear. Without loss of generality, assume 
that pj lies on the line segment between pi and pk- This means that Xj lies on 
the c-segment with respect to xq between xt and Xk- Thus by O Proposition 
5.1] and the compactness of f2 , 

max{-c(x,Xi) + c(x , x t ), -c(x, x k ) +c(x ,x k )} 

> -c(x,Xj) + c(x ,Xj) + C d M (x l ,x k ) 2 d M (x,x ) 2 - jd M (x,x ) 3 

for some C > 0, 7 > and all d^/ (x,xq) sufficiently small. Since Xi 7^ xi, if 
dM (x, xq) < r for r > sufficiently small we will have C" d M (xi, Xk) dM (x, xq) — 
jdM (x,xq) > 0, and hence either — c(x,Xi) + c(xo,Xi) 7^ —c(x,Xj) +c(xo,Xj) 
or — c(x, Xk) + c(xo,Xk) 7^ — c(x,Xj) + c(xo,Xj) as desired. □ 

The following lemma is the most computationally intensive of the paper. In 
calculating the difference quotient of G % , wc will encounter intersections of sets 
of the form W^ A \ Wj or Wj \ W^ A . What we show is that an intersection 
of these sets with two different indices has [i measure that decays like o(|Ad|) 
as I Ad I — s> 0. This will allow us to eliminate most of the terms in the final 
expression of -§^-G % . 
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Now, under assumption QMTW+D , each of the sets W Ad \ W 3 or Wj \ W Ad 
are essentially the differences of two nested, strongly convex sets which decrease 
to (d — 1) -dimensional sets as | Ad| — > 0. If we consider the intersection of two 
such "generalized annuli" and the limiting sets intersect transversally, the decay 
rate of o(|Ad|) is readily imagined. In the case when the limiting sets intersect 
tangentially, we can apply Lemma 15.31 above to obtain the desired decay The 
reader who is satisfied with the preceeding explanation may wish to skip over 
the proof of the following lemma. 

Lemma 5.4. If c, fi, tt, and fi satisfy all of the conditions in Section^ then 
for any 1 < j ^ k < K , 

fi {{W Ad \ Wj) n (W L Ad \ W k )) = o(|Ad|), 
ft {{Wj \ w Ad ) n {w k \ W k Ad )) = (|Ad|). 

as |Ad| — > 0. Here, the rate of decay o(|Ad|) may depend on d £ M. K . 

Proof. By conditions (|Twist[) , ( |Nondeg I , the boundedness of Vt, and since ft 
is absolutely continuous with respect to dVoU/ with bounded density, for any 
measurable set E C fi we have 

n(E)<\[E] £i \ Cn , 

where is the volume measure induced on T£.M by the Riemannian metric 
on M. Hence, it will be sufficient to prove 

(\W A % \ \w 3 ] s ) n ([w Ad ]^ \ \w k ] s ) \ cn = o(|Ad|) 

(the second equality in the statement of the lemma follows by a symmetric 
argument). Now since Ad 1 = 0, this implies W Ad = Wi, hence we can assume 
that j, k i. We then define 

Xi(p) : = exp%.{p) 

C{,»(p) : = c(Xi(p),xi) -c(Xi(p),Xi) 

for any 1 < I < K, so that 





= {pe [n] g( 


-Cj,i(p) < logd 3 - 


logd 1 } 






= {pe [«]*, 


-Cfc,i(p) < lo S dfc - 


logcf} 




[^%, 




-Cj,i{p) < log (cP - 


hA#) - 


logcf} 






-c k Ap)<^g{d k 


f Ad k ) - 


-logcf} 



(5.2) 

If either Ad? or Ad k are non-positive, we would have 

{[w Ad ] s \[w 3 ] s )n([w Ad ] s \[w k ] 
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hence we may assume both quantities are strictly positive. 

We will now show that the desired rate of decay holds on a neighborhood 
of each point contained in [cWj]=. H [9Wfc]=. More precisely, we show that for 
each p £ [dWj]~ fl [cWfc] 2 . , there exists a neighborhood M p depending on d and 
the cost function c such that 

= o(|Ad|) (5.3) 

as | Ad| — > 0. There are a number of different cases to work through, depending 
on how the two sets [9Wj] 2 . and [9Wfe] 2 . intersect. It is understood that when- 
ever we take a new "small neighborhood" of a point, it will be contained in any 
previous such neighborhoods. Also, in each of the following cases, we assume 
by a translation that 

o g [dWi\ Si n [dW k ] Si . 

Case 1: (Intersection points in the interior of [f2] 2 .) 

Suppose that is in the interior of [f2] s .. We have two subcases. 
Case la: (Tangential intersection) 

First assume that the normal vectors to [9Wj] 2 . and [<9Wfc] 2 . at are parallel. 
Since is in the interior of in some small neighborhood of 0, 

[dWj] Si C{ P £ [0] 2j | _S,.<(p) = logd j - logd*} , 
[dW k ] Si C{pe [Sl] Si | -c M (p) = logd' - logd*} 

and —D5j t i and —Dc k .i (nonzero by (| Twist [I ) are normal vectors to [9Hj-] s . and 
[9Wfc] 2 . respectively. Thus by Lemma IS~3l we can find some radius r > such 
that 

Wx i nW ii n5 r (0) = {0}. (5.4) 

By translating and rotating coordinates (denoted by (p 1 , . . . ,p d )), we may as- 
sume that the (d — l)-dimcnsional hypcrplanc normal to — D5j j(0) is given by 
p d = 0, with —Dcj t i(0) in the positive p d direction. Since c is diffcrentiablc, 
[cWj] 2 is given as the graph of a real valued C 1 function of p' = (p 1 , . . . ,p d ~ 1 ) 
on a small neighborhood of 0, call this function p(p'). Consider the transforma- 
tion given by 

p i y (p',p d - p(p')), 

which maps [cWj] 2 onto the hyperplane P d := {p d = 0} in this neighborhood. 
It is easy to see that the differential of this map is the identity matrix at the ori- 
gin, and so it has nondegencrate Jacobian determinant on a small neighborhood 
of the origin. Thus, we may assume that 

[dW^ c P d , 
[Wi-L c { P d < 0} (5.5) 



([< d ] 2i \ [WjQ n ([^ Ad ] 2i \ [W k ]- X ) nAfp 
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in some small neighborhood of while (|5.4p continues to hold. 

We will now bound the height of ([W J Ad ]_ \ [Wj] £ ) n ([W^ Ad ] _ \ [W k ] s 

in the p d direction in a small neighborhood of the origin, in terms of |Ad|. 
Recalling that 

(-Dcj ti (p),e d ) = \-Dc jt i(p)\ > 

for p G [<9Wj]- by (j Twist [) . we can find some sufficiently small neighborhood 
B ro (0) of the origin on which 

(-Dcj,i(p),e d ) > ^(-Dc M (p'),e d ). 

Now hx any p G ([V^ Ad ] _ \ [Wj] s ^ H B ro (0) (so in particular, p' G [dWj]-.), 
then since Ad- 7 > and by (|5.5j) we see that 

(p - p ', e d ) > 0. 

Then we calculate, 
Ad j 

—— + o(|Ad|) = log(cP + Ad?) - logd» 
dP 

»1 

{-Dc j:i {tp+{l-t)p'),(p-p'))dt 



\p-p'\ [ (-Dc j , i (tp+Q.-t)j/),e d )dt 
Jo 



> \{-D~c 3 Ap')^d)\p-p'\ 
= \\-Dc hl {p')\\p-p'\ 
>\p-p'\ 

by (|Twist[) and the compactness of ft. In particular, for some C" > depending 
on d and c, 

|p-p'|<C"|Ad| 



for | Ad | small enough, and any p G f[ W/ d ] -. \ [Wj] s .j H B ro (0). Thus 

(V Ad ] £i \ [W 3 ] s ) nB ro (0) c {|/| < C"|Ad|} nB ro (0). (5.6) 
We now claim that for jVq := B ro (0), 

diam (V(([Vlf % \ [W^) H \ \W k ] Sf )) n M.) = o(l) (5.7) 
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as |Ad| — > 0, where n d (-) is orthogonal projection P d . Suppose not, then there 
is a constant Kq > 0, a sequence of Ad„ with |Ad„| — > as n — > oo, and a 
sequence of points 



Pl,n, P2,n G 



W Ad„ 



such that \p'i n — P2, n \ > Ko for all n. By the boundedness of A/"o, we may pass 
to subsequences and assume that pi tTi — > p\ and p2, n — > P2 for some points p\ 
and pi as n — > oo. By the above calculations, we see that 

\pi,n-p[, n \ <<?"|Ad„| 
\ P 2,n-p' 2 J<C"\Ad n \ 

for each n. Hence, by letting n — > oo we see that p\ = p[ and p 2 = p'21 an d thus 

\pi — P2I > ^o- 

However, by the continuity of c, it is clear that both p\ and p 2 must be contained 
in [dWj] m . ("1 [9Wfc] 2 . flJVo which contradicts 

Finally, by §M the set ([Hf d ]_ ( \ D ([^ Ad ]^ \ [W fc ] s J nA^o is 

contained in a cylinder with (c? — 1) -dimensional base 

* d (([ W f d ] Sl \ Wh) n ([^^ \ [W*lx 4 )) nA/ "° 

and height C"|Ad|, and by (|5.7[) we see this set must have measure o(|Ad|) as 
|Ad| —J- 0, thus obtaining 

Case lb: (Nontangential intersection) 

Now suppose the normal vectors to [cWj] 2 . and [9Wfc]-. at are not parallel. 
Again, translating, rotating, and straightening out [<?Wj] 2 . near 0, we assume 
that —Dcj t i(0) is in the positive p d direction, and both of the inclusions (|5.5[) 
hold in a small neighborhood of the origin. 

Since the intersection point is in the interior of [fi] 2 . , by ( |Reg[ ) and (|Twist[) , 
there is a unique (d— l)-dimensional tangent hyperplane T [dWk] s ., to [9Wfc] s . 
at and it does not coincide with P d . Since they do not coincide, the intersection 
of To [cWfc] 2 . and P d form an (d — 2) -dimensional linear subspace. Pick a basis 
{vi, . . . , Vd~ 2 } for this subspace, then add two vectors v and w to the collection 
such that {vi, . . . , Vd-2, v } is a basis for P d and {v\, . . . , Vd-2,w} is a basis for 
To [cWfc] 2 . . Clearly, {v\, . . . , v^-2, v, w} is a linearly independent collection of 
d vectors. If we define a linear transformation L by 

L{vi)=v u \<l<d-2, 
L(v) = v, 

L(u;) = (0,...,0,l), 

this L is invertible, maps P d onto itself, and maps To [cWfc] 2 . onto an (d — 1) 
dimensional linear subspace orthogonal to P d . By an additional rotation, we 
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may assume this subspace is P d_1 := =0}. Since [<9Wfc] s . is locally 

the graph of a C 1 function on P d_1 , we may straighten [c?Wfc] 5 . on a small 
neighborhood Wo of so it coincides with P d_1 (as in the first part of the proof 
of Case la), and this transformation will map P d to itself (in particular, we 
will still have n Mo C P d ). Thus, we may assume that in some small 

neighborhood of the origin, (|5.5j) holds and [cWfc] 2 . = P d_1 , while [W-fc] 2 . C 
{p 1 < 0}. Additionally, by restricting to a smaller neighborhood if necessary 
(which we still write as Mo), we obtain by similar reasoning to Case la above 
that 

(K Ad ] 2l \ m s ) n ([W£*\ ai \ [W*]. 4 ) nM, 
c {I/- 1 ) < C"|Ad|, \ P d \ < C"\Ad\} nMo, 

which has measure of order o(|Ad|), hence we obtain (|5.3|) again. 

Case 2: (Intersection points on [<9f2] 2 .) 

Suppose now that € [<9f2] 2 .. Again we consider two subcases. 
Case 2a: ("Fake boundary") 

Suppose that —5^(0) ^ logd- 7 — logd* (thus —5,^(0) < logd- 7 — logd 1 by the 
definition of Wj). By the continuity of —Cj,%, we can find a small neighborhood 
Mo of the origin that depends on c and d on which both —Cjj < log (dp + Ad- 7 ) — 
logd' and —Cji < \ogdP — logcf whenever |Ad| is sufficiently small, i.e. 

(K Ad ] Si \[^-] Ss )nA^ = 0, 

and we immediately see (|5.3j) . 

A symmetric argument holds if — 6^(0) ^ \ogd k — logd 4 . 
Case 2b: (Tangential intersection) 

Now suppose that —5^(0) = logd- 7 — logd* = \ogd k — logd' = — 6/^,(0), and 
the normal vectors —Dcj^(0) and —Dck,i(0) are parallel. Then, note that since 
[0]-. is convex and bounded by assumption, there exists a Lipschitz mapping 
that straightens the boundary of [Wj] £ . H [f2] 2 . in a small neighborhood of the 
origin, with a controlled Lipschitz norm depending only on f2, c, and d. As 
such, we may apply the same proof as Case la above to obtain (|5.3[) (note that 
Lemma 15.31 is still valid when the point xo E dfl). 

Case 2c: (Nontangential intersection) 

For the final case, suppose that — Cj»(0) = logd- 7 — logd' = logd fc — logd 4 = 
— Cfc,i(0), and the normal vectors —Dcj,i(0) and —Dck,i(0) are not parallel. Since 
c is C 1 up to the boundary of fl, we may extend both — £.,•,, and — Ck.i outside 
of to a small neighborhood of in a C 1 manner. Then we may apply the 
same proof as in Case lb to obtain (|5.3p (note that we do not require the use 
of Lemma 15.31 in this case) . 

Now we can show the desired global equality. For each point p <G [9Wj]-, H 
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[<9Wfc] 2 ., there exists a neighborhood M p that corresponds to p via one of the 
above steps, and we can use the compactness of [9Wj-] s . PI [<9Wfc] 2 . to extract a 
finite cover M Pl , . . . M PN . It is easy to see that each of the sets [W^* 1 ] _ \ [Wj]- . 
and [W^ 1 ] _ \ [W/s] 2 . are contained in some neighborhoods of [cWj]-, and 
[9W)t] s . respectively, whose diameter decreases to zero with |Ad|. Thus we 
obtain 

\(\wf % \ [wjQ n ([wn^ \ m] s ,)\ cn 

N 

1=1 
= o(\Ad\) 

as | Ad | — > 0, as desired. □ 

With this lemma in hand, we are finally ready to show the differentiability 
of G 1 , and establish an upper bound for the magnitude of the derivative. The 
calculations have been split into two separate propositions, due to length. 

Proposition 5.5. Assume c, fi, 0, and Cl satisfy all of the conditions in Sec- 
tion [3J Then, 



G*(d + Ad) - G'(d) = J2 ^ A f dk + 

k=l 



\k,Ad : / I(Xi(p))\-Dc k ,i(p)\ |" X |det DX^p)!^" 1 (p) 

•* £^fc,d,Ad 



whe 



and (using some of the same notation as the proof of Lemma \5.4 
Xi{p) : = exp%.{p), 
Ck,i{p) ■= c(Xi{p),Xk) - c{X l (p),x i ), 



£4,d,Ad 



' [ P e [ii] Si | -c M (p) = \o g d k - logd*} n [n^ fe 
{ P G [Q] Si | -Cfc,i(p) - log + Ad fe ) - logcf } n [n, #fc Wi 



Ad k > 
Ad fc < 0. 



Proof. Wc can calculate, in the notation of Lemma [4 



G t (d + Ad)-G i (d) = / i f]W ; 



Ad 



\/=l / 



/^(Vd+AdC^)) - H(Vd(Xi)) 



A* (Vd+Ad,i) - A 4 



a* fn< d \ n - a* f n ^ \ fW d ) 

\i=l k=l ) \k=l 1=1 ) 



18 



Now note that 

/ K K 



Mn< d \nH 

\i=i fe=i / 

= J\J (f)wf d \w k 



\k=l \l=l 

K / K 



K 



K 



fc=i \;=i 



\i=i 



E 



f| iy Ad \w k An ( f) iy Ad \ w fe2 ] n 
a=i / \i=i / 



A" 



n 



For any fixed index k, we have that f]f =1 W t Ad \ W k C V^ Ad \ W fc , hence by 
Lemma T5. 41 above, all terms except the first one in the last line of calculations 
above has a rate of decay o(|Ad|) as |Ad| — » 0. Since we obtain a similar 
expression for 



K 



K 



» n w * \ n w i 

\k=i i=i 



Ad 



we find 

G l (d + Ad) -G*(d) 
/' 



E 

fe=i 

AT 

E 

fc=i 



K 



oflAdl 



f| w/ Ad n (v^ Ad \ w k ) \ /i f) w l n \ 



Ad\ 



o(| Ad|) 



Now, fix an index 1 < k < K. First assume that Ad k > 0, hence W fe \W fc Ad 
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Then, by making the change of variable p = —Dc(x,Xi), we see that 



K 



K 



fx p| wf* n (wt A \ W k ) - m f| Wi n \ w, 



Ad\ 

k ) 



K 



f) w t Ad n (wf d \ w k ) 

A=ik 



n^ k w^n(w^\w k ) 



[r\^ k w^n(w^\w k )} 



I(x)dVol M {x) 



IiXiip^detDX^ldp 



= (*)• 

Recalling f|5 . 2[) . we see that 



A' 

n < d n (^ Ad \ w k) 

l^k 



log d k - log d* < -Ck ti (p) < log + Ad*)- log d* 



Then since —Ck.iip) is a Lipschitz function on [f2] s ., by applying the coarea 
formula, we obtain 



(*) 



log (d J +Ad fe )-logd' 
log d fc —log d 1 



I{X t (p) ) | -D~c K% (p) | - 1 |dct DX t (p) IdH"" 1 (p) 



{-c fc , i ( P )=t}n[n i7!fc ^ 

/(X, (p))| -Dcfc.i (p) I" 1 1 dct DX t (p) | dU d - 1 (p) ] (log (d k + Ad k ) - log d k ) + o{ | Ad|) 



£4., 



1 



J(X i (p))|-D5 fc , l (p)|- i |detDX i (p)|^ d - 1 (p) 1 Ad fc +o(|Ad|) 



^fc,d,Ad 

k 



(|Ad|) 



(5.8) 



_ /3»,fc,AdA d fc 
d fc 

as |Ad| -> 0. 

If AcZ fe < we have Wj^ d \Wfc = instead, and following similar calculations 
we again obtain the expression (|5.8[) . By summing over 1 < k < X, the proof 
is completed. □ 



We are now ready to provide the necessary upper bound on -Jj- 
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Proposition 5.6. G % is differentiable in d l and 



,d KCsu Pn l[H d ^(d[n} St )] 

l dd~ 1 - 

for each 1 < i < K , where 

|det (-D5c(x,Xi))| 



C := max sup ■ 



i<k#<K aea \(-DDc(x,xi)) l (-Dc(x,xi)+Dc(x,x k ))\ Taiii 

Proof. Fix d > 0. Now it is easy to sec that for any A > we have G l (d) = 
G l (\d). We have used here (jTwistl) . ( |NondegD , and the fact that I is bounded 
above to move the o(t) terms outside of the summation. Thus, applying Propo- 
sition [53] above, we calculate for any t with \t\ < d l , 

G i (d\...,d i + t,...,d K )-G l (d) 

d 1 d % d K d i 

= G 4 (-^, ...../' - Gi ( d ) 

y d l +t d l +t' x ' 

= G\d x - f£t + o{t), ...,d i ,...,d K -^ r t + o(t)) - G*(d) 



Eft,fc,Ad t / d k \ 
— ]k- +°(*) 



fc=i v 
t * 

= -^E^d t +o(t) (5.9) 
fe=i 

as f — > 0, and with Ad t := (— -pt, . . . , 0, . . . , ——t) with the in the ith compo- 
nent. Then, by applying the dominated convergence theorem we first see that 
^j-G l exists. 

In order to obtain the claimed inequality, first note that by the convexity of 
[fi] 2 . (by assumption) and Theorem 14.41 the level sets of — c are the boundaries 

of convex subsets of [Q]^. . Thus by Lemma I5T21 we have 

U d - X (f/ fe ,d,Ad t ) < U d ~ x (d[£l] Si ) 

for any 1 < k < K and \t\ small. Additionally, by the definition of Xi and C, 
and (|Twist[) and (Nondcgl, we can calculate for any 1 < k < K, 

l-Dcfc^r^detDXiO?)! < C, 
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hence by combining with (|5.9p we obtain 



d_ & = . G i {d 1 , . . . ,d l + t, . . . , d K ) - G'(d) 

A' 



fe=l 

as claimed. □ 

6 An upper bound on the number of steps 

Recall that when we initialize the scheme, we require an initial choice of d° such 
that G 4 (d°) < fi + 5 for all 2 < i < K. This can be attained if we make the 
choice of 

( 4 = 1 (6.1) 

where 

M:= max sup e" c(a; ' 2fc)+c(a; ' 2l) + 1. 

In fact, with this choice we have G l (d ) = for all 2 < i < K. Now notice that 
in each step, we either do not change the value of each d l or decrease it, so we 
may assume these values as upper bounds for our choices of d l throughout the 
algorithm. Now, recalling that the value of d 1 does not change (hence d 1 = 1 
throughout), and that at every step our choice of </>d remains in the set 4>^, we 
can use (|3.4[) to calculate 

< h ~ K8 
<fi-(K- 1)6 

K 



i=2 
A 



<l-J2G l (d) 

i=1 

= G 1 (d). 



Now if we write 



A:= min inf e -c(*,*i,)+c(»,xi) 

Kk<K i6f! 



we can see from definition that Va(xi) = if d l < A for any i ^ 1, and by 
Lemma H31 we would have G : (d) = 0. Hence there is a lower bound of 

d l > A (6.2) 
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for any d 1 with i 1. 

Now, say we are at the (n, i 



l)-st step of the algorithm, and we must 



decrease d l n 4-1 to In this case we can combine Proposition l5.6l with (|6.2j) 

to obtain, 

< G l (0„, 4 ) - G i {(f) nt i-i) 

= G l (d 



n.i— 1 ) * • * J u n,i- 1 ) 



dd, G 



< sup 

ATC sup n / 



< 



A 



> ^n,i— l) ^ '(^n,i— 1' • • • ' ^"n,i— 1> * 
(<^n,i-l — ^ri,i-l) 



i a n,i-l) 



in other words, we obtain a strictly positive lower bound for the magnitude of 
the decrease in the parameter d l . 

However, since only strictly positive values of d l are admissible and the 
parameters d l can only decrease, we see that each d % can only be decreased a 
finite number of times. The worst case scenario is when only one parameter is 
updated per iteration of the scheme, hence the maximum number of iterations 
is given by K times the maximum number of times each d l can be updated. 
Thus by the initial choice (|6.ip . we see that an upper bound on the number of 
iterations n e that this algorithm can take is given by the claimed bound (|3 . 3|) : 



n F < K 



KCM su Po / 
SA 



max W 1 - 1 (d[Q] 9 ) 

Ki<K v x%> 
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